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QUASI- THREE-DIMENSIONAL FLOW SOLUTION BY 
MERIDIONAL PLANE ANALYSIS 
by Theodore Katsanis and William D. McNally 
Lewis Research Center 

ABSTRACT 

A computer program has been developed to obtain subsonic or shock- 
free transonic, nonviscous flow analysis on the hub- shroud mid- channel 
flow surface of a turbomachine. The analysis may be for any annular 
passage, with or without blades. The blades may be fixed or rotating and 
may be twisted and leaned. The flow may be axial, radial or mixed. 

Blade surface velocities over the entire blade are approximated based 
on the rate of change of angular momentum. This gives a 3-D flow picture 
based on a 2-D analysis. 

The paper discusses the method used for the program and shows ex- 
amples of the type of passages and blade rows which can be analyzed. Also, 
some numerical examples are given to show how the program can be used 
for practical assistance in design of blading, annular passages, and 
annular diffusers. 

f 

INTRODUCTION 

The design of blades for compressors and turbines ideally requires 
analysis methods which include unsteady, rotational, three-dimensional, 
and viscous effects. Clearly, such solutions are impractical at the present 
time, even on the largest and fastest computers. The usual approach at 
present is to analyze only steady flows and to separate inviscid solutions 
from viscous solutions. At present., inviscid analyses usually involve a 
combination of two-dimensional solutions oh intersecting families of stream 
surfaces to obtain what is called a quasi-three-dimensional solution. 
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Since there are several choices of two-dimensional surfaces to analyze, 
and many ways of combining them, there are many approaches to obtaining 
a quasi-three-dimensional solution. Most two-dimensional solutions are 
either on a blade-to-blade surface of revolution (Wu's surface, (1)'% see 
fig. 1) or on the mid-channel stream surface between two blades (Wu’s Sg 
surface). Another type of two-dimensional solution can be obtained on a 
passage cross-sectional surface (normal to the flow). 

In this paper a solution to the equations of flow on the mid- channel Sg 
surface is presented. This solution surface is chosen when the turbo- 
machine under consideration has significant variation in flow properties in 
the hub-shroud direction. A solution on the mid-channel surface will show 
this variation. The solution can be obtained either by the quasi- orthogonal 
method, which solves the velocity- gradient equation from hub to shroud on 
the mid- channel flow surface (2), or by a finite-difference method, which 
solves a finite-difference equation for stream function on the same flow 
surface. The quasi-orthogonal method is efficient in many cases and can 
obtain solutions into the transonic regime. However, there is a convergence 
problem when aspect ratios are above 1 or when the passage has high wall 
curvature compared to the passage hei^t. For such cases, the most 
promising method is the finite-difference solution, but this solution is 
limited to completely subsonic flows. 

The method described in this paper uses the finite-difference method, 
followed by velocity gradient method, if necessary. When there is tran- 
sonic flow, the finite difference method is used to obtain a reduced mass 
flow subsonic solution for flow angles and streamline curvatures. This 

Numbers in parentheses designate references at end of paper. 
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information is then used with the velocity gradient method to obtain the 
transonic solution. 

A computer program, MERIDL, has been written to perform these cal- 
culations (3 and 4). This program is written for axial-, radial-, or mixed- 
flow in any annular passage, with or without blades. Upstream and down^ 
stream flow conditions can vary from hub to shroud. The solution is for 
compressible, shock-free flow, or incompressible flow. Provision is made 
for an approximate correction for loss of stagnation pressure through the 
blade row. The blade row may be either fixed or rotating and may be twisted 
and leaned. The blades can have high aspect ratio and arbitrary thickness 
distribution. 

The solution obtained by this program also provides the information 
necessary for a more detailed flow analysis on blade-to-blade surfaces 
(fig. 1). A useful program for this purpose is TSONIC (5). Information 
needed to prepare all the input for TSONIC is calculated and printed by 
MERIDL, 

To give an idea of the type of problem that can be anlyzed with MERIDL, 
a few illustrations are given. The first illustration is an axial compressor 
blade row with high aspect ratio and high twist, shown in figure 2. 

Another type of problem is shown in figure 3. The entire annular 
passage here could be analyzed as 4 cases. First, is the inlet mixed flow 
bladeless passage, second the axial stator, third the axial rotor, and finally 
the mixed flow diffuser. There can be whirl in the ducts, and the whirl, 
temperature and stagnation pressure can vary from hub to shroud. Each 
of the 4 cases would be analyzed separately, and each case would have to 
be matched to the adjacent case. 

Figure 4 shows another row type where this method applies. This is 
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a radial turbine. The stator and rotor would be analyzed separately. Since 
this stator would be expected to have uniform flow from hub to shroud, it could 
best be analyzed by a blade-to-blade program such as TSONlC (5). The rotor 
could be analyzed with MERIDL. The rotor could have either radial of swept 
blades^ MERIDL has not been coded to analyze a case with splitter blades. 

This paper describes the MERIDL program computer requirements and 
the computer input and output. The metiiod of analysis is explained briefly, 
with a list of the basic assumptions. Finally, three numerical examples 
are given to show the results that can be obtained. 

COMPUTER PROGRAM 

The computer program MERIDL has been implemented on the NASA 
Lewis time sharing IBM-TSS/360-67 computer as well as the NASA Lewis 
Univac 1106 computer. Storage for program variables requires 60,000 
words for a 41-axial by 21-radial grid; storage for the program code is 
18,000 words additional. Run times range from 3 to 15 minutes on the IBM 
360-67 equipment, depending upon the mesh size used and the compress- 
ibility of the flow. The MERIDL program is available through COSMIC at 
the University of Georgia, and is documented thoroughly in (3 and 4). 

Input to the program consists primarily of blade section and hub- shroud 
geometry as well as upstream andcbwnstream flow-variable boundary con- 
ditions. The major input is summarized on figure 5. Hub and shroud 
geometries are both given by z vs. r spline coordinates. Blade geometry 
is described by a series of blade sections from hub to shroud on blade- 
to-blade surfaces of revolution. Meanline ( 0 ) and blade thickness coor- 
dinates are given as functions of z and r to define these blade sections. 
Upstream and downstream flow boundary conditions are given at a number 
of locations from hub to shroud. Upstream of the blade row, total tern- 
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perature, total pressure, and inlet whirl (tangential velocity) are required, 
while downstream total pressure (or loss) and outlet whirl are given. These 
conditions, along with mass flow, define the inviscid problem to be solved. 

Output from the program can be obtained at any points within the solu- 
tion region at which it is requested by the user. Ordinarily output is re- 
quested along user-designated streamlines, the locations of which are cal- 
culated by the program. It can also be obtained on the finite-difference 
solution mesh. Output is calculated after each major iteration of the solu- 
tion and after any velocity gradient solutions. The user can control the 
number of iterations at which he desires it to be printed. Figure 6 shows 
the most important output which usually requested. 

There are many types of output which can be given, but the most 
common variables printed at a requested set of points are the following: 


Streamline Coordinates . . 
Meridional Velocity .... 
Tangential Velocity .... 

Relative Velocity 

Meridional Flow Angle . , 
Blade-to-Blade Flow Angle 
Critical Velocity Ratio . , 
Streamline Curvature 
Blade Surface Velocities . 




• • • 

w, 

w 

or 

^ 

• w/w„ 

Wp (within blade regions) 


See figure 7 and the list of Symbols for the definition of these variables. 
Extensive plots can also be obtained of both the geometrical and fluid input 
variables and of the major output variables. The principal outputs plotted 
are the calculated streamline pattera in the meridional plane as well as 
meridional and blade surface velocities along these streamlines. 

When a reasonable flow pattern has been achieved by MERIDL on the 
mid-channel flow surface, more detailed and accurate blade surface veloc- 
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ities can be obtained on blade-to-blade surfaces which follow the direction of 
the calculated meridional streamlines. Several blade elements or stream 
channels can be analyzed from hub to shroud using the TSONIC program (5) 
to build a quasi-three-dimensional solution for the full flow passage. MERIDL 
calculates and prints most of the input required for the TSONIC blade-to- 
blade program along each of the calculated stream channels (fig. 6). 

METHOD OF ANALYSIS 

BASIC ASSUMPTIONS - It is desired to determine the flow distribution 
through an annular duct, or through a stationary or rotating cascade of blades 
on a mid-channel hub-shroud stream surface. The following simplifying 
assumptions are used in deriving the equations and in obtaining a solution; 

(1) The flow relative to a blade is steady. 

(2) The flow is axi-symmetric where there is no blade. 

(3) The fluid is a perfect gas with constant specific heat Cp. 

(4) The fluid is a nonviscous gas. 

(5) There is no heat transfer. 

(6) The mid- channel surface is a stream surface which has the same 
shape as the blade mean camber surface, except near the leading and trailing 
edges, where an arbitrary correction is made to match the free- stream flow. 

(7) The only forces are those due to momentum and pressure gradient. 

(8) The velocity varies linearly between blade surfaces. 

(9) The relative stagnation pressure loss is known through the blade row. 

The flow may be axial, mixed, or radial. There may be a variation of 

whirl, stagnation pressure, and stagnation temperature from hub to shroud, 
both upstream and downstream of the blade row. The blade row m^y be 
either fixed or rotating, with leaned and twisted blades. Within the given 
assumptions, no terms are omitted from the equations. 
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SOLUTION BY COMBINATION OF METHODS - A flow analysis on the 
meridional flow surface can be obtained either by the velocity-gradient 
method or by the finite-difference method. The finite-difference solution of 
the stream function equation is limited to subsonic flow, whereas the velocity- 
gradient method by itself is limited to relatively low-aspect- ratio blades. 

Of these two methods the most accurate solution is obtained by the finite- 
difference technique, so that this method is used where possible (i. e, , for 
subsonic flow). In cases which have locally supersonic flow, the finite- 
difference solution is first obtained at a reduced mass flow for which the 
flow field is completely subsonic. The streamline curvatures and flow angles 
throughout the passage which: are obtained from this solution provide the 
information necessary to obtain an approximate velocity- gradient solution 
at full mass flow, even with high aspect ratio and high wall curvature. 

SUBSONIC STREAM-FUNCTION SOLUTION - The stream-function equa- 
tion is a partial differential equation on a mid- channel hub- shroud stream 
surface (see assumption 6). This equation is in one unknown (the stream 


function u) as a function of two variables, rrand z (see fig. 7). 


3^u . 3^u 1 3u /l 


0r- 
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F 

^ dr p do 



The derivatives of the stream function satisfy the equations 
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( 2 ) 

Equation (1) is derived from Wu’s equation (eq. (107a), (1)) for the stream 
function on what he calls an S 2 surface. Equation (1) is nonlinear but can be 
solved iteratively by the finite-difference method when the flow is completely 
subsonic. 

A finite region (as indicated in fig. 8) is considered for the solution of 
equation (1). It is assumed that the upstream and downstream boundaries 
are sufficiently far from the blade so as to have a negligible effect on the 
solution. Equation (1) is elliptic for subsonic flow. Therefore, when the 
flow is entirely subsonic, equation (1) can be solved when proper boundary 
conditions are specified on the entire boundary of the region. These con- 
ditions are the values of the stream function on all four boundaries. The 
stream function has the value 0 on the hub and 1 at the shroud. The value 
of the stream function on the upstream and downstream boundaries can be 
calculated if the stagnation pressure, stagnation temperature, and whirl dis- 
tribution from hub to shroud are specified upstream and downstream of the 
blade. 

The numerical solution of equation (I) is obtained by the fiiiite-difference 
method. The finite-difference grid is an orthogonal mesh which is gener^ 
ated by the program, using the method reported in reference (6). Figure 9 
illustrates the orthogonal finite-difference mesh. 

The finite-difference equations are nonlinear since the original equa- 
tion (1) is nonlinear. These equations can be solved iteratively. On the 
first iteration an initial density is assumed; this linearizes some of the 
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terms. The remaining nonlinear terms are omitted for the first iteration 
so that the finite-difference equations are entirely linearized. These 
linearized equations are then solved to obtain the first approximate solu- 
tion for stream function. This solution provides information used to obtain 
a better estimate of the density and an estimate of the other nonlinear terms. 
The equations are the solved again to obtain an improved solution. This 
process is repeated, and by iteration a final converged solution can be ob- 
tained if the flow is subsonic. 

For each step of this iteration, the linearized finite-difference equations 
must be solved. The method used to solve the equations is successive over- 
relaxation (8) with an optimum overrelaxation factor. Since this is also 
an iterative method, we have two levels of iteration. 

After the stream function is obtained, the velocity distribution is ob- 
tained by numerical partial differentiation of the stream function and by 
using equation (2). The details of the numerical procedure and programming 
technique are described in reference (4). 

TRANSONIC VELOCITY- GRADIENT APPROXIMATE SOLUTION - For 
the case where there is locally supersonic flow, equation (1) is no longer 
elliptic in the entire region but is hyperbolic in the region of supersonic 
flow (7). This changes the boundary conditions and means that there will 
probably be shock losses in going from supersonic to subsonic flow. The 
finite-difference method cannot be used with locally supersonic flow. 
However, an approximate solution can be obtained by getting a reduced- 
flow solution with the finite-difference method and extending this to the full 
flow by using the velocity- gradient method. This technique is described in 
reference (5). 
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The velocity-gradient equation is 

dW = (aW + b + — + d cos dt + — + Wf ■ (3) 

V W / W 

The coefficients, a, b, c, and d, are given by different expressions in the 
blade region, in the upstream region, and in the downstream region. These 
coefficients are given as follows: 

Blade region coefficients: 

^ ^ cos^^ cos (a _ sin^£cosjP ^ ^ |3 cos ^ ^ 


c = 0 • 

d = 0 

Upstream- region coefficients: 

a - cos(a - (?) 
b = 0 

c = - k-J£.A p cos(a- <p) + cos 

V p2 /[ ^c 

d = — sinCa - (p) 

dm 



/dW. ^ 

b = cos/3 — — sin(a - 2oJsinPcos<^ + rcos jSi— + 2co sin a 

dm Vdm / 




11 


Downstream- region coefficients: 


a - cos(of - <p) 


b = 0 


c = - 


w 


- ior 


0 


Mo ■ 


ct)r 


cos(a - (p) + 


TVn\ + O^T' 
O 


> (6) 


cos cpl 


dW 


d = 


m 


dm 


sin(a - (p) 




Finally, in all three regions, we have 


e = C T! - o> dX - C dT" +^3^dp" 
pi P v^t1 


P’ 




f = - Rdp^^ 

2C p” 2T" 




(7) 


Equations (3) to (7) are derived in reference (3). Equation (3) is solved 
as an initial-value problem, where the initial value of W is specified at the 
hub for any given mesh line running from hub to tip. By finding several solu- 
tions for varying values of W at the hub, a solution satisfying continuity will 
be found; that is, the solution will satisfy 


y^pW rB cos(a - <p) cos j3 dt = w 


( 8 ) 


When equation (3) has been solved, subject to satisfying equation (8), 
for every hub-to-tip mesh line in the region, the entire velocity distribution 
is obtained. 

BLADE SURFACE VELOCITIES - The solution which is obtained by 
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either the finite-difference or velocity- gradient method is for the mid-channel 
surface between the blades. Of greafef interest are the blade surface veloc- 
ities from which boundary layer growth and losses can be calculated. Blade 
surface velocities can be estimated from the mid- channel solution since the 
blade loading is dependent on the rate of change of whirl. Blade surface ve- 
locities can be calculated by assuming a linear variation of velocity between 
blade surfaces. In reference (3) the following equation for calculating blade 
surface velocities is derived: 


B 

^ 2 dm 


B d(rV.) 

^tr ^ 

mia 2 dm 


(9) 


NUMERICAL EXAMPLES 

AXIAL-FLOW COMPRESSOR ROTOR - Flow was analyzed on the mid- 
channel flow surface of an axial-flow rotor designed with the computer pro- 
gram of (9). The design pressure ratio is 1. 275; the inlet hub-tip radius 
ratio, 0.5; the aspect ratio, 1.5; the tip solidity, 1,0; and the tip relative 
Mach number at the inlet, 0. 9. Although tip relative Mach number is near 
sonic, there are no locally supersonic regions on the meridional flow sur- 
face. This permitted the calculation and comparison of both the finite dif- 
ference and velocity gradient solutions. 

Figures 10, 11, and 12 show plots from MERIDL for this example. 
Figure 10 shows the input hub, shroud, and blade leading and trailing edges 
in the meridional plane. Design blade sections are also indicated. There 
were 861 mesh points: 41 in the axial direction and 21 in the radial di- 
rection. The input hub, mean and tip blade sections are shown in figure 11. 
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Calculated mid- channel and surface velocities obtained from MERIDL 
are shown in figure 12 for the hub, 50 percent- and tip streamlines. 

A comparison, for the mid-channel velocities only, between the stream- 
function solution and the more approximate velocity gradient solution for 
the same weight flow is. also shown on figure 12. The dashed-line velocity 
gradient solutions deviate very little from the stream function solution, and 
only do so in the blade regions where maximum turning is taking place. 

More details on this example, including both input and output, are given 
in (3). 

AXIAL -TO- RADIAL DIFFUSER DUCT - This duct is the diffuser for an 
axial flow turbine. Because of space limitation the flow was turned to the 
radial direction in a short distance. The flow angle ^ varies from 15° at 
the hub to 10° at the shroud, at the turbine rotor exit. The maximum Mach 
number is less than 0.4, so that the flow is fully subsonic. 

Figure 13 shows the inner and outer walls of the diffuser, with the in- 
put spline points marked. There were 400 mesh points; 25 along hub and 
shroud and 16 across the passage from hub to shroud, as shown in figure 14. 
The streamlines obtained by the MERIDL finite-difference solution for incre- 
ments of 10 percent flow are shown in figure 15, The velocity distribution 
at the inner wall, 50 percent streamline, and outer wall are shown in fig- 
ure 16. The static pressure distribution along the inner and outer walls, 
as a function of the true streamline distance (not meridional distance) are 
shown in figure 17. These can be used in a boundary layer analysis to 
locate possible separation points. 

FLARED AXIAL STATOR - This example is a stator for one stage of. a 
multi-stage axial turbine. The hub, shroud, and blade leading and trailing 
edges are shown in figure 18. The blade sections at hub, mean, and tip 
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are shown in figure 19. 

There is constant whirl entering the blade row. Because of the curva- 
ture and change of radius of the hub and shroud, the whirl angle changes 
rapidly just ahead of the blade leading edge. Small changes in hub or shroud 
change these flow angles drastically, which can lead to large incidence losses. 
The MERIDL finite-difference solution was used to calculate the flow angle 
just ahead of leading edge for several hub and shroud contours. The flow 
angles for the final configuration (fig. 18) are compared with the blade angles 
in figure 20. 

SYMBOLS 

a coefficient in velocity- gradient equation (3) 

B tangential space between blades, rad 

b coefficient in velocity- gradient equation (?) 

1 

Cp specific heat at constant pressure, J/ (kg) (K) 
c coefficient in velocity- gradient equation (3) 
d coefficient in velocity- gradient equation (3) 
e coefficient in velocity-gradient equation (3) 

F vector normal to mid- channel stream surface and proportional to 

tangential pressure gradient, N/1^ 

f coefficient in velocity- gradient equation (3) 

2 2 

1 rothalpy, C^TI - ojX, meters /sec 

m meridional streamline distance, meters 

2 

p pressure, N/meters 

R gas constant, J/(kg) (K) 

r radius from axis of rotation, meters 

r radius of curvature of meridional streamline, meters 
c 

T temperature, K 
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t distance along orthogonal mesh lines in direction across flow 
u normalized stream function 
W fluid velocity relative to blade^ meters/sec 
w mass flow, kg/sec 
z axial coordinate, meters 

a angle betweem meridional streamline and axis of rotation, rad; see 
figure 7 

angle between relative velocity vector and meridional plane, rad; see 
figure 7 

? coefficient in stream-function equation, defined in equation (1) 

& relative angular coordinate, rad; see figure 7 

2 

X prerotation, (rV^)^, meters /sec 

^ coefficient in stream- function equation, defined in equation (1) 

3 

p density, kg/meter 

(p angle between s- coordinate line and axis of rotation, rad; see figure 4 
oj rotational speed, rad/sec; see figure 7 
Subscripts: 
cr critical 
i inlet 

1 blade surface facing direction of positive rotation 
m component in direction of meridional streamline 
mid mid- channel blade to blade 
o outlet 

r component in radial direction 
tr blade surface facing direction of negative rotation 
component in tangential direction 


0 



Superscripts: 

’ absolute stagnation condition 
relative stagnation condition 
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Figure 2, - Axial con^ressor 





CD-9967-Ol 

Figure 3. - Axial turbine with mixed flow annular passage. 



Figure 4 -RadiaJ inflow turtine. 




Figures. - MERIDL program input. 
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^STREAMLINE COORDINATES z, r 
MERIDIONAL FLOW ANGLE a 
STREAMLINE CURVATURE 
VELOCITY COMPONENTS Wq 
RELATIVE VELOCITY W 
BLADE-TO-BLADE FLOW ANGLE, p 
BLADE SURFACE VELOCITIES Wj, 
CRITICAL VELOCITY RATIO W/W^.^ 


Figure 6. - MERIDL program output. 






Figure 7. - Cylindrical coordinate system and velocity components. 


SHROUD-* ^SOLUTION REGION 



Figures. -Mid-channel surface solution region. 
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Figure 9. - Orthogonal mesh on the solution region. 



Z DIRECTION 

Figure 10. * Axial flow compressor rotor with design blade sections. 



Figure 
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(cl SHROUD STREAMUNE. 


Figure 12. - Meridional and blade surface velocities for com- 
pressor rotor. 




Figure 14. - Solution orthogonal mesh for diffuser duct. 
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COORDINATE, FT 
ic) OUTER WALL. 

Figure 16. - Velocity distribution 
for diffuser duct. 
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(a) INNER WALL STREAMLINE. 



STREAMLINE DISTANCE FROM UPSTREAM BOUNDARY. FT 



(b) OUTER WALL STREAMLINE. 

Figure 17. - Static pressure distribution for diffuser duct 
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<a) HUB. 



(c)TIP. 

Figure 19. - Axial stator blade sections. 
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Figure 20. - Inlet flow angle compared with blade angle for 
axial flow stator. 



